Load data

# RNA SEQ DATA
load("/Users/senosam/Documents/Massion_lab/RNASeq_summary/rnaseq.RData")
ls_preprocessed <- preprocess_rna(path_rnaseq = '/Users/senosam/Documents/Massion_lab/RNASeq_summary/rnaseq.RData', correct_batch = T, correct_gender = T)
p_all <- p_all[which(p_all$Vantage_ID %in% colnames(ls_preprocessed$vsd_mat)),]
p_all <- na.omit(p_all[match(k, p_all$pt_ID),])
rna_all <- ls_preprocessed$rna_all[,c(1:7,na.omit(match(p_all$Vantage_ID, colnames(ls_preprocessed$rna_all)[8:ncol(ls_preprocessed$rna_all)])))]

pData_rnaseq <- ls_preprocessed$pData_rnaseq[na.omit(match(p_all$pt_ID, ls_preprocessed$pData_rnaseq$pt_ID)),]
counts_all <- ls_preprocessed$counts_all[,na.omit(match(p_all$Vantage_ID, colnames(ls_preprocessed$counts_all)))]
vsd_mat <- ls_preprocessed$vsd_mat[,na.omit(match(p_all$Vantage_ID, colnames(ls_preprocessed$vsd_mat)))]
ls_preprocessed2 <- list(p_all=na.omit(p_all), 
                         rna_all=rna_all, 
                         pData_rnaseq=pData_rnaseq, 
                         counts_all=counts_all, 
                         vsd_mat=vsd_mat)

p_all_cl <- ls_preprocessed2$p_all
cl2 <- data.frame(cl[match(p_all$pt_ID, cl$pt_ID),])
cl2$pt_ID <- as.character(cl2$pt_ID)
p_all_cl <- p_all_cl %>% 
              inner_join(., cl2, "pt_ID") %>%
              mutate_all(as.character)
p_all_cl['Gender'] <- ls_preprocessed2$pData_rnaseq$Gender

Compare CyTOF cell types

x <- data.frame(CyTOF_prcnt)
x['TP_53'] <- as.factor(cl$TP_53)
x <- na.omit(x)
x <- reshape2::melt(x,  id.vars = c('TP_53'))

violin_deconv <- function(x, ct) {
  x_sub <- x[which(x$variable == ct),]
  dp <- ggplot(x_sub, aes(x=TP_53, y=value, fill=TP_53)) + 
    geom_violin(trim=FALSE) +
    geom_boxplot(width=0.1, fill="white") +
    labs(title=paste0(ct, " - CyTOF"),x="pt group", y = "cell %") +
    scale_fill_brewer(palette="Dark2") + theme_minimal() +
    ggsignif::geom_signif(comparisons = list(c("0", "1")),   
                map_signif_level=TRUE)
  dp
}

xcell_ct <- colnames(CyTOF_prcnt)
for (i in xcell_ct){
  plot(violin_deconv(x, ct = i))
}

## Warning in wilcox.test.default(c(0.01404193449648, 0.0582096474953618,
## 0.000117145400871562, : cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0067722681359045, 0.0702690166975881, :
## cannot compute exact p-value with ties

Compare xCell cell types

x <- data.frame(RNA_xcell)
x['TP_53'] <- as.factor(cl$TP_53)
x <- na.omit(x)
x <- reshape2::melt(x,  id.vars = c('TP_53'))

violin_deconv <- function(x, ct) {
  x_sub <- x[which(x$variable == ct),]
  dp <- ggplot(x_sub, aes(x=TP_53, y=value, fill=TP_53)) + 
    geom_violin(trim=FALSE) +
    geom_boxplot(width=0.1, fill="white") +
    labs(title=paste0(ct, " - XCell deconvolution"),x="pt group", y = "Enrichment score") +
    scale_fill_brewer(palette="Dark2") + theme_minimal() +
    ggsignif::geom_signif(comparisons = list(c("0", "1")),   
                map_signif_level=TRUE)
  dp
}

xcell_ct <- colnames(RNA_xcell)
for (i in xcell_ct){
  plot(violin_deconv(x, ct = i))
}
## Warning in wilcox.test.default(c(0, 0.0085, 0, 0, 0, 0, 0, 0, 0, 0, 0.0062, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0396, 0, 0, 0.2611, 0.0438, 0, 0.0482, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.2829, 0.0549, 0, 0, 0, 0.0068, 0, 0.0776, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0224, 0, 0, 0.0312, 0, 0, 0.0464, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0029, 0.0153, 0.0539, 0, 0.0223, 6e-04, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0668, 0.0823, 0, 0.1576, 4e-04, 0.0798, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0118, 0.0307, 0, 0, 0.0485, 0.0106, 0.0014, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.1122, 0.004, 0, 0.1344, 0.0127, 0, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0774, 0.0299, 0, 0.1108, 0.0061, 0, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0451, 0.0147, 0, 0.0884, 0, 0, 0.0364, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0578, 0.0148, 0, 0.0357, 0, 0, 0, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0054, 3e-04, 0, 0.0034, 0, 0.0013, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0, 0.0036, 0.0078, 0.0083, 0.0115, 0.0421, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0081, 0, 0.012, 0, 0, 0, 0, 0, 0, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0.1264, 0.1048, 0, 0.0905, 0.0655, 0.0453, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0306, 0.0219, 0.0346, 0, 0.1207, 0.0651, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0128, 0.0114, 0.004, 0, 0, 0.0012, 0, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.1926, 0.1741, 0.0551, 0.2679, 0.0486, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0043, 0.0046, 1e-04, 0, 0, 0, 0, 0, 0, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0541, 0.091, 0.2541, 0.1425, 0.1426, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0.2287, 0.3249, 0.0653, 0.2771, 0.1111, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0014, 0, 0.033, 0, 0.0085, 0, 0, 0, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.3187, 0.2078, 0.1312, 0.4575, 0.0889, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0018, 0, 0, 0.007, 0, 0, 0, 0, 0, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0019, 5e-04, 0.019, 0.0083, 0.0077, 0.0105, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0, 0.0506, 0.0196, 0.0242, 0.0063, 0.0716, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.1354, 0.0114, 0, 0.0238, 0, 0.1263, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0081, 0.1566, 0.0244, 0.036, 0.0119, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0625, 0.017, 0.013, 0.0081, 0.0093, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0619, 0.1094, 0.0208, 0.064, 0, 0.0426, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.057, 0.052, 0.0127, 0.0441, 0.0501, 0.0444, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0066, 0.0054, 0, 0.0108, 1e-04, 0.0026, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0047, 0.0081, 0, 0.0135, 0, 0.0087, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0, 0, 0, 0.0831, 0, 0, 0, 0, 0, 0.002, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0127, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.0472, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0234, 0, 0, 0, 0, 0, 0.054, 0, 0, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0204, 0.0122, 0.0144, 0.0184, 0.0032, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.028, 0.0456, 0, 0.0263, 0, 0.007, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0026, 0.0012, 7e-04, 0.0029, 0.0039, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.0024, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0, 0.1068, 0, 0.0166, 0, 0.1044, 0.0115, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0306, 0, 0, 0.0488, 0, 0.0866, 0, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0011, 0, 0.0353, 0, 0.0725, 0.0062, 0.0458, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0041, 0.0071, 0, 0.0034, 0.007, 0, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.133, 0.0177, 0.0705, 0.0498, 0.0636, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 4e-04, 0.0123, 0, 0.0088, 0.0036, 0.0059, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.002, 0, 0, 0, 0, 0, 0, 0, 0.0217, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0, 0.0641, 0, 0.046, 0, 0.0393, 0.0542, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.054, 0.0177, 0.0077, 0.0162, 0.0102, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 5e-04, 0, 0, 0, 0, 0, 0.0538, 0, 0, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0337, 0, 0, 0.0379, 0, 0, 0.0093, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0477, 0.076, 0.0098, 0.0498, 0.0176, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0.5017, 0.4452, 0.2353, 0.1457, 0.0375, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.1447, 0.0936, 0.0221, 0.1745, 0.014, : cannot
## compute exact p-value with ties

## Warning in wilcox.test.default(c(0.0808, 0.1142, 0.037, 0.0985, 0.027, 0.0334, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0084, 0, 0, 0.018, 0.0064, 0, 0.0081, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0.0019, 0, 0, 0, 0, 0, 0.0158, 0, 0.0059, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0, 0, 0, 0, 0.0127, 0, 0, 1e-04, 0, 0.0037, :
## cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(0.2106, 0.2538, 0.0602, 0.2725, 0.0798, :
## cannot compute exact p-value with ties

Compare SILA score

x_sila <- data.frame(cbind('TP_53' = cl$TP_53, 'SILA' = CDE_TMA36$SILA))
x_sila$TP_53 <- as.factor(x_sila$TP_53)
x_sila <- reshape2::melt(x_sila,  id.vars = c('TP_53'))
x_sila$TP_53 <- as.factor(x_sila$TP_53)

ggplot(x_sila, aes(x=TP_53, y=value, fill=TP_53)) + 
    geom_violin(trim=FALSE) +
    geom_boxplot(width=0.1, fill="white") +
    labs(title=paste0("SILA score"),x="pt group", y = "SILA score") +
    scale_fill_brewer(palette="Dark2") + theme_minimal() +
    ggsignif::geom_signif(comparisons = list(c("0", "1")),   
                map_signif_level=TRUE)
## Warning in wilcox.test.default(c(0.61, 0.695, 0.479, 0.665, 0.761, 0.417, :
## cannot compute exact p-value with ties

Compare HLA-DR expression

In ECC (CyTOF)

x <- data.frame(cbind('TP_53' = cl$TP_53, 'HLA-DR' = CyTOF_exp$Epi_174Yb_HLA.DR))
x$TP_53 <- as.factor(x$TP_53)
x <- reshape2::melt(x,  id.vars = c('TP_53'))
x$TP_53 <- as.factor(x$TP_53)

ggplot(x, aes(x=TP_53, y=value, fill=TP_53)) + 
    geom_violin(trim=FALSE) +
    geom_boxplot(width=0.1, fill="white") +
    labs(title=paste0("HLA-DR expression in ECC"),x="pt group", y = "HLA-DR expression") +
    scale_fill_brewer(palette="Dark2") + theme_minimal() +
    ggsignif::geom_signif(comparisons = list(c("0", "1")),   
                map_signif_level=TRUE)
## Warning: Removed 12 rows containing non-finite values (stat_ydensity).
## Warning: Removed 12 rows containing non-finite values (stat_boxplot).
## Warning: Removed 12 rows containing non-finite values (stat_signif).

In Myeloid cells (CyTOF)

x <- data.frame(cbind('TP_53' = cl$TP_53, 'HLA-DR' = CyTOF_exp$Mye_174Yb_HLA.DR))
x$TP_53 <- as.factor(x$TP_53)
x <- reshape2::melt(x,  id.vars = c('TP_53'))
x$TP_53 <- as.factor(x$TP_53)

ggplot(x, aes(x=TP_53, y=value, fill=TP_53)) + 
    geom_violin(trim=FALSE) +
    geom_boxplot(width=0.1, fill="white") +
    labs(title=paste0("HLA-DR expression in Myeloid cells"),x="pt group", y = "HLA-DR expression") +
    scale_fill_brewer(palette="Dark2") + theme_minimal() +
    ggsignif::geom_signif(comparisons = list(c("0", "1")),   
                map_signif_level=TRUE)
## Warning: Removed 12 rows containing non-finite values (stat_ydensity).
## Warning: Removed 12 rows containing non-finite values (stat_boxplot).
## Warning: Removed 12 rows containing non-finite values (stat_signif).

In Mesenchymal cells (CyTOF)

x <- data.frame(cbind('TP_53' = cl$TP_53, 'HLA-DR' = CyTOF_exp$FMes_174Yb_HLA.DR))
x$TP_53 <- as.factor(x$TP_53)
x <- reshape2::melt(x,  id.vars = c('TP_53'))
x$TP_53 <- as.factor(x$TP_53)

ggplot(x, aes(x=TP_53, y=value, fill=TP_53)) + 
    geom_violin(trim=FALSE) +
    geom_boxplot(width=0.1, fill="white") +
    labs(title=paste0("HLA-DR expression in FMes"),x="pt group", y = "HLA-DR expression") +
    scale_fill_brewer(palette="Dark2") + theme_minimal() +
    ggsignif::geom_signif(comparisons = list(c("0", "1")),   
                map_signif_level=TRUE)
## Warning: Removed 12 rows containing non-finite values (stat_ydensity).
## Warning: Removed 12 rows containing non-finite values (stat_boxplot).
## Warning: Removed 12 rows containing non-finite values (stat_signif).

In RNA Seq (main genes)

hladr_a <- 'ENSG00000204287'
hladr_b1 <- 'ENSG00000196126'
col11a <- 'ENSG00000060718'

HLA_mut <- cbind('TP_53' = cl$TP_53, 
               'HLADR_A' = RNA_top12K_E[,grep(hladr_a, colnames(RNA_top12K_E))], 
               'HLADR_B1' = RNA_top12K_E[,grep(hladr_b1, colnames(RNA_top12K_E))],
               'COL11A1' = RNA_top12K_E[,grep(col11a, colnames(RNA_top12K_E))])
HLA_mut <- data.frame(na.omit(HLA_mut))
HLA_mut <- reshape2::melt(HLA_mut,  id.vars = c('TP_53'))
HLA_mut$TP_53 <- as.factor(HLA_mut$TP_53) 
HLA_mut$value <- as.numeric(HLA_mut$value)

violin_clust <- function(x, ct) {
  x_sub <- x[which(x$variable == ct),]
  dp <- ggplot(x_sub, aes(x=TP_53, y=value, fill=TP_53)) + 
    geom_violin(trim=FALSE) +
    geom_boxplot(width=0.1, fill="white") +
    labs(title=paste0(ct, " gene expression"),x="", y = "VST") +
    scale_fill_brewer(palette="Dark2") + theme_minimal() +
    ggsignif::geom_signif(comparisons = list(c("0", "1")),   
                map_signif_level=TRUE)
  dp
}

gene <- unique(HLA_mut$variable)
for (i in gene){
  plot(violin_clust(HLA_mut, ct = i))
}
## Warning in wilcox.test.default(c(4.42573767657876, 13.7344965547443,
## 13.8165079801241, : cannot compute exact p-value with ties

## Warning in wilcox.test.default(c(4.09877639859111, 14.3175852872763,
## 14.9297465989775, : cannot compute exact p-value with ties

DEG

DE_res_1v2 <- DE_analysis(ls_preprocessed2, 
           GeneBased=FALSE, 
           pDataBased=FALSE,
           NewCondition=TRUE,
           NewCondition_df = p_all_cl,
           cond_nm='TP_53',
           two_levels=c('0','1'),
           reference = '0',
           correct_gender=TRUE)
## Unlist done
## # A tibble: 52 x 5
##    Vantage_ID  pt_ID Batch TP_53 Gender
##    <chr>       <chr> <chr> <chr> <fct> 
##  1 R3388_YZ_46 11424 1     1     Male  
##  2 R3388_YZ_1  11601 1     0     Male  
##  3 R3388_YZ_2  11646 1     0     Male  
##  4 R3388_YZ_4  11652 1     0     Male  
##  5 R3388_YZ_44 11817 1     1     Female
##  6 R3388_YZ_3  11820 1     0     Male  
##  7 R4163_YZ_28 11840 2     1     Male  
##  8 R3388_YZ_5  11855 1     0     Male  
##  9 R3388_YZ_43 11938 1     0     Male  
## 10 R3388_YZ_59 11957 1     0     Male  
## # … with 42 more rows
## # A tibble: 52 x 5
##    Vantage_ID  pt_ID Batch Condition Gender
##    <chr>       <chr> <chr> <chr>     <fct> 
##  1 R3388_YZ_46 11424 1     1         Male  
##  2 R3388_YZ_1  11601 1     0         Male  
##  3 R3388_YZ_2  11646 1     0         Male  
##  4 R3388_YZ_4  11652 1     0         Male  
##  5 R3388_YZ_44 11817 1     1         Female
##  6 R3388_YZ_3  11820 1     0         Male  
##  7 R4163_YZ_28 11840 2     1         Male  
##  8 R3388_YZ_5  11855 1     0         Male  
##  9 R3388_YZ_43 11938 1     0         Male  
## 10 R3388_YZ_59 11957 1     0         Male  
## # … with 42 more rows
## # A tibble: 52 x 5
##    Vantage_ID  pt_ID Batch Condition Gender
##    <chr>       <chr> <chr> <chr>     <fct> 
##  1 R3388_YZ_46 11424 1     1         Male  
##  2 R3388_YZ_1  11601 1     0         Male  
##  3 R3388_YZ_2  11646 1     0         Male  
##  4 R3388_YZ_4  11652 1     0         Male  
##  5 R3388_YZ_44 11817 1     1         Female
##  6 R3388_YZ_3  11820 1     0         Male  
##  7 R4163_YZ_28 11840 2     1         Male  
##  8 R3388_YZ_5  11855 1     0         Male  
##  9 R3388_YZ_43 11938 1     0         Male  
## 10 R3388_YZ_59 11957 1     0         Male  
## # … with 42 more rows
## # A tibble: 52 x 5
##    Vantage_ID  pt_ID Batch Condition Gender
##    <chr>       <chr> <chr> <chr>     <fct> 
##  1 R3388_YZ_46 11424 1     1         Male  
##  2 R3388_YZ_1  11601 1     0         Male  
##  3 R3388_YZ_2  11646 1     0         Male  
##  4 R3388_YZ_4  11652 1     0         Male  
##  5 R3388_YZ_44 11817 1     1         Female
##  6 R3388_YZ_3  11820 1     0         Male  
##  7 R4163_YZ_28 11840 2     1         Male  
##  8 R3388_YZ_5  11855 1     0         Male  
##  9 R3388_YZ_43 11938 1     0         Male  
## 10 R3388_YZ_59 11957 1     0         Male  
## # … with 42 more rows
## Labeling done
## Filtering done
## Warning in DESeqDataSet(se, design = design, ignoreRank): some variables in
## design formula are characters, converting to factors
## Design done
## Warning in DESeqDataSet(se, design = design, ignoreRank): some variables in
## design formula are characters, converting to factors
## Warning: Setting row names on a tibble is deprecated.
## vsd symbols done
## using pre-existing size factors
## estimating dispersions
## gene-wise dispersion estimates
## mean-dispersion relationship
## final dispersion estimates
## fitting model and testing
## -- replacing outliers and refitting for 2217 genes
## -- DESeq argument 'minReplicatesForReplace' = 7 
## -- original counts are preserved in counts(dds)
## estimating dispersions
## fitting model and testing
## DESeq done
## using 'normal' for LFC shrinkage, the Normal prior from Love et al (2014).
## 
## Note that type='apeglm' and type='ashr' have shown to have less bias than type='normal'.
## See ?lfcShrink for more details on shrinkage type, and the DESeq2 vignette.
## Reference: https://doi.org/10.1093/bioinformatics/bty895
## res symbols done
## list done
heatmap_200(DE_res_1v2$res_df, DE_res_1v2$vsd_mat_sym, DE_res_1v2$meta_data, DE_res_1v2$pData_rnaseq)

volcano_plot(DE_res_1v2$res_df, gene=NULL, p_title='TP53: WT vs Mut')

fgsea_res_1v2 <- fgsea_analysis(DE_res_1v2)
## `summarise()` ungrouping output (override with `.groups` argument)
## Warning in fgsea(pathways = gmtPathways(pthw_path), stats = ranks, nperm
## = 1000): You are trying to run fgseaSimple. It is recommended to use
## fgseaMultilevel. To run fgseaMultilevel, you need to remove the nperm argument
## in the fgsea function call.
## Warning in preparePathwaysAndStats(pathways, stats, minSize, maxSize, gseaParam, : There are ties in the preranked stats (0.01% of the list).
## The order of those tied genes will be arbitrary, which may produce unexpected results.
## Warning in fgsea(pathways = gmtPathways(pthw_path), stats = ranks, nperm
## = 1000): You are trying to run fgseaSimple. It is recommended to use
## fgseaMultilevel. To run fgseaMultilevel, you need to remove the nperm argument
## in the fgsea function call.
## Warning in preparePathwaysAndStats(pathways, stats, minSize, maxSize, gseaParam, : There are ties in the preranked stats (0.01% of the list).
## The order of those tied genes will be arbitrary, which may produce unexpected results.
## Warning in fgsea(pathways = gmtPathways(pthw_path), stats = ranks, nperm
## = 1000): You are trying to run fgseaSimple. It is recommended to use
## fgseaMultilevel. To run fgseaMultilevel, you need to remove the nperm argument
## in the fgsea function call.
## Warning in preparePathwaysAndStats(pathways, stats, minSize, maxSize, gseaParam, : There are ties in the preranked stats (0.01% of the list).
## The order of those tied genes will be arbitrary, which may produce unexpected results.
## Warning in fgsea(pathways = gmtPathways(pthw_path), stats = ranks, nperm
## = 1000): You are trying to run fgseaSimple. It is recommended to use
## fgseaMultilevel. To run fgseaMultilevel, you need to remove the nperm argument
## in the fgsea function call.
## Warning in preparePathwaysAndStats(pathways, stats, minSize, maxSize, gseaParam, : There are ties in the preranked stats (0.01% of the list).
## The order of those tied genes will be arbitrary, which may produce unexpected results.
## Warning in fgsea(pathways = gmtPathways(pthw_path), stats = ranks, nperm
## = 1000): You are trying to run fgseaSimple. It is recommended to use
## fgseaMultilevel. To run fgseaMultilevel, you need to remove the nperm argument
## in the fgsea function call.
## Warning in preparePathwaysAndStats(pathways, stats, minSize, maxSize, gseaParam, : There are ties in the preranked stats (0.01% of the list).
## The order of those tied genes will be arbitrary, which may produce unexpected results.
## Warning in fgsea(pathways = gmtPathways(pthw_path), stats = ranks, nperm
## = 1000): You are trying to run fgseaSimple. It is recommended to use
## fgseaMultilevel. To run fgseaMultilevel, you need to remove the nperm argument
## in the fgsea function call.
## Warning in preparePathwaysAndStats(pathways, stats, minSize, maxSize, gseaParam, : There are ties in the preranked stats (0.01% of the list).
## The order of those tied genes will be arbitrary, which may produce unexpected results.
## Warning in fgsea(pathways = gmtPathways(pthw_path), stats = ranks, nperm
## = 1000): You are trying to run fgseaSimple. It is recommended to use
## fgseaMultilevel. To run fgseaMultilevel, you need to remove the nperm argument
## in the fgsea function call.
## Warning in preparePathwaysAndStats(pathways, stats, minSize, maxSize, gseaParam, : There are ties in the preranked stats (0.01% of the list).
## The order of those tied genes will be arbitrary, which may produce unexpected results.
## Warning in fgsea(pathways = gmtPathways(pthw_path), stats = ranks, nperm
## = 1000): You are trying to run fgseaSimple. It is recommended to use
## fgseaMultilevel. To run fgseaMultilevel, you need to remove the nperm argument
## in the fgsea function call.
## Warning in preparePathwaysAndStats(pathways, stats, minSize, maxSize, gseaParam, : There are ties in the preranked stats (0.01% of the list).
## The order of those tied genes will be arbitrary, which may produce unexpected results.
## Warning in fgsea(pathways = gmtPathways(pthw_path), stats = ranks, nperm
## = 1000): You are trying to run fgseaSimple. It is recommended to use
## fgseaMultilevel. To run fgseaMultilevel, you need to remove the nperm argument
## in the fgsea function call.
## Warning in preparePathwaysAndStats(pathways, stats, minSize, maxSize, gseaParam, : There are ties in the preranked stats (0.01% of the list).
## The order of those tied genes will be arbitrary, which may produce unexpected results.
fgsea_res <- fgsea_res_1v2
cond_nm <- 'TP53: WT vs Mut'
fgsea_plot(fgsea_res$res_hm, pathways_title='Hallmark', condition_name= cond_nm)

## # A tibble: 30 x 8
##    pathway                    pval    padj     ES   NES nMoreExtreme  size state
##    <chr>                     <dbl>   <dbl>  <dbl> <dbl>        <dbl> <int> <chr>
##  1 HALLMARK_G2M_CHECKPOINT 0.00174 0.00542 -0.677 -3.28            0   187 down 
##  2 HALLMARK_E2F_TARGETS    0.00175 0.00542 -0.650 -3.16            0   191 down 
##  3 HALLMARK_EPITHELIAL_ME… 0.00176 0.00542 -0.554 -2.68            0   189 down 
##  4 HALLMARK_MTORC1_SIGNAL… 0.00178 0.00542 -0.493 -2.39            0   190 down 
##  5 HALLMARK_ALLOGRAFT_REJ… 0.00174 0.00542 -0.491 -2.38            0   187 down 
##  6 HALLMARK_MITOTIC_SPIND… 0.00175 0.00542 -0.483 -2.35            0   195 down 
##  7 HALLMARK_INFLAMMATORY_… 0.00176 0.00542 -0.455 -2.21            0   193 down 
##  8 HALLMARK_INTERFERON_GA… 0.00175 0.00542 -0.443 -2.16            0   194 down 
##  9 HALLMARK_IL6_JAK_STAT3… 0.00179 0.00542 -0.506 -2.16            0    82 down 
## 10 HALLMARK_TNFA_SIGNALIN… 0.00175 0.00542 -0.425 -2.07            0   195 down 
## # … with 20 more rows
fgsea_plot(fgsea_res$res_c1, pathways_title='C1 positional genes', condition_name= cond_nm)

## # A tibble: 30 x 8
##    pathway     pval   padj     ES   NES nMoreExtreme  size state
##    <chr>      <dbl>  <dbl>  <dbl> <dbl>        <dbl> <int> <chr>
##  1 chr19q13 0.0016  0.0103 -0.493 -2.73            0   769 down 
##  2 chr16q22 0.00230 0.0103  0.572  2.66            0   132 up   
##  3 chr4q35  0.00220 0.0103  0.625  2.38            0    43 up   
##  4 chr5q13  0.00220 0.0103  0.566  2.37            0    70 up   
##  5 MT       0.00220 0.0103  0.672  2.34            0    31 up   
##  6 chr16q12 0.00224 0.0103  0.600  2.33            0    48 up   
##  7 chr3q21  0.00175 0.0103 -0.523 -2.25            0    90 down 
##  8 chr2q21  0.00182 0.0103 -0.563 -2.23            0    60 down 
##  9 chr5q21  0.00212 0.0103  0.657  2.20            0    26 up   
## 10 chr1p33  0.00223 0.0103  0.588  2.18            0    40 up   
## # … with 20 more rows
fgsea_plot(fgsea_res$res_c2, pathways_title='C2 curated genes', condition_name= cond_nm)

## # A tibble: 30 x 8
##    pathway                     pval   padj     ES   NES nMoreExtreme  size state
##    <chr>                      <dbl>  <dbl>  <dbl> <dbl>        <dbl> <int> <chr>
##  1 SOTIRIOU_BREAST_CANCER_… 0.00180 0.0135 -0.818 -3.79            0   138 down 
##  2 ROSTY_CERVICAL_CANCER_P… 0.00181 0.0135 -0.820 -3.78            0   130 down 
##  3 SHEDDEN_LUNG_CANCER_POO… 0.00173 0.0135 -0.679 -3.58            0   426 down 
##  4 KONG_E2F3_TARGETS        0.00187 0.0135 -0.810 -3.49            0    88 down 
##  5 DUTERTRE_ESTRADIOL_RESP… 0.00180 0.0135 -0.680 -3.48            0   305 down 
##  6 WHITEFORD_PEDIATRIC_CAN… 0.00189 0.0135 -0.785 -3.47            0   108 down 
##  7 CROONQUIST_IL6_DEPRIVAT… 0.00184 0.0135 -0.791 -3.44            0    91 down 
##  8 ZHOU_CELL_CYCLE_GENES_I… 0.00187 0.0135 -0.766 -3.44            0   116 down 
##  9 FLORIO_NEOCORTEX_BASAL_… 0.00180 0.0135 -0.713 -3.43            0   171 down 
## 10 KOBAYASHI_EGFR_SIGNALIN… 0.00183 0.0135 -0.680 -3.41            0   239 down 
## # … with 20 more rows
fgsea_plot(fgsea_res$res_c3, pathways_title='C3 regulatory target genes', condition_name= cond_nm)

## # A tibble: 30 x 8
##    pathway                 pval   padj     ES   NES nMoreExtreme  size state
##    <chr>                  <dbl>  <dbl>  <dbl> <dbl>        <dbl> <int> <chr>
##  1 HSD17B8_TARGET_GENES 0.00159 0.0350 -0.572 -3.07            0   539 down 
##  2 E2F_Q6               0.00177 0.0350 -0.509 -2.48            0   212 down 
##  3 E2F1_Q6              0.00176 0.0350 -0.504 -2.46            0   215 down 
##  4 E2F4DP1_01           0.00176 0.0350 -0.499 -2.45            0   219 down 
##  5 E2F_Q4               0.00176 0.0350 -0.500 -2.44            0   215 down 
##  6 SGCGSSAAA_E2F1DP2_01 0.00176 0.0350 -0.523 -2.44            0   153 down 
##  7 E2F_Q4_01            0.00180 0.0350 -0.496 -2.42            0   216 down 
##  8 E2F1DP1RB_01         0.00176 0.0350 -0.495 -2.41            0   213 down 
##  9 E2F_02               0.00179 0.0350 -0.486 -2.37            0   218 down 
## 10 E2F1DP1_01           0.00179 0.0350 -0.485 -2.37            0   218 down 
## # … with 20 more rows
fgsea_plot(fgsea_res$res_c4, pathways_title='C4 cancer', condition_name= cond_nm)

## # A tibble: 30 x 8
##    pathway       pval    padj     ES   NES nMoreExtreme  size state
##    <chr>        <dbl>   <dbl>  <dbl> <dbl>        <dbl> <int> <chr>
##  1 GNF2_CCNA2 0.00198 0.00852 -0.867 -3.53            0    65 down 
##  2 MODULE_54  0.00174 0.00852 -0.697 -3.48            0   241 down 
##  3 GNF2_PCNA  0.00198 0.00852 -0.854 -3.47            0    65 down 
##  4 GNF2_CCNB2 0.00190 0.00852 -0.886 -3.47            0    54 down 
##  5 GNF2_CDC20 0.00187 0.00852 -0.877 -3.43            0    53 down 
##  6 GNF2_CDC2  0.00192 0.00852 -0.852 -3.42            0    59 down 
##  7 GNF2_MCM4  0.00187 0.00852 -0.872 -3.40            0    52 down 
##  8 GNF2_BUB1B 0.00189 0.00852 -0.870 -3.36            0    49 down 
##  9 GNF2_CENPF 0.00189 0.00852 -0.836 -3.35            0    58 down 
## 10 GNF2_HMMR  0.00195 0.00852 -0.881 -3.33            0    46 down 
## # … with 20 more rows
fgsea_plot(fgsea_res$res_c5, pathways_title='C5 GO genes', condition_name= cond_nm)

## # A tibble: 30 x 8
##    pathway                     pval   padj     ES   NES nMoreExtreme  size state
##    <chr>                      <dbl>  <dbl>  <dbl> <dbl>        <dbl> <int> <chr>
##  1 GO_SISTER_CHROMATID_SEG… 0.00175 0.0304 -0.588 -2.78            0   171 down 
##  2 GO_MITOTIC_SISTER_CHROM… 0.00180 0.0304 -0.593 -2.72            0   142 down 
##  3 GO_DNA_DEPENDENT_DNA_RE… 0.00178 0.0304 -0.586 -2.69            0   138 down 
##  4 GO_NUCLEAR_CHROMOSOME_S… 0.00174 0.0304 -0.548 -2.68            0   230 down 
##  5 GO_CHROMOSOME_CENTROMER… 0.00175 0.0304 -0.563 -2.66            0   172 down 
##  6 GO_REGULATION_OF_CHROMO… 0.00185 0.0304 -0.624 -2.66            0    92 down 
##  7 GO_CHROMOSOME_SEGREGATI… 0.00169 0.0304 -0.527 -2.66            0   287 down 
##  8 GO_COLLAGEN_FIBRIL_ORGA… 0.00189 0.0304 -0.695 -2.64            0    51 down 
##  9 GO_CELL_CYCLE_DNA_REPLI… 0.00190 0.0304 -0.668 -2.64            0    61 down 
## 10 GO_CHROMOSOMAL_REGION    0.00169 0.0304 -0.514 -2.61            0   298 down 
## # … with 20 more rows
fgsea_plot(fgsea_res$res_c6, pathways_title='C6 oncogenic', condition_name= cond_nm)

## # A tibble: 30 x 8
##    pathway                   pval   padj     ES   NES nMoreExtreme  size state
##    <chr>                    <dbl>  <dbl>  <dbl> <dbl>        <dbl> <int> <chr>
##  1 RB_P107_DN.V1_UP       0.00180 0.0130 -0.597 -2.75            0   128 down 
##  2 CSR_LATE_UP.V1_UP      0.00181 0.0130 -0.531 -2.50            0   156 down 
##  3 PRC2_EED_UP.V1_DN      0.00176 0.0130 -0.483 -2.30            0   179 down 
##  4 BMI1_DN_MEL18_DN.V1_UP 0.00181 0.0130 -0.497 -2.29            0   135 down 
##  5 PRC2_EZH2_UP.V1_DN     0.00179 0.0130 -0.481 -2.28            0   173 down 
##  6 RPS14_DN.V1_DN         0.00180 0.0130 -0.480 -2.28            0   172 down 
##  7 LEF1_UP.V1_UP          0.00176 0.0130 -0.453 -2.16            0   179 down 
##  8 SNF5_DN.V1_UP          0.00182 0.0130 -0.459 -2.15            0   152 down 
##  9 ATF2_UP.V1_UP          0.00181 0.0130 -0.446 -2.10            0   165 down 
## 10 E2F1_UP.V1_UP          0.00179 0.0130 -0.440 -2.09            0   173 down 
## # … with 20 more rows
fgsea_plot(fgsea_res$res_c7, pathways_title='C7 immunologic', condition_name= cond_nm)

## # A tibble: 30 x 8
##    pathway                    pval    padj     ES   NES nMoreExtreme  size state
##    <chr>                     <dbl>   <dbl>  <dbl> <dbl>        <dbl> <int> <chr>
##  1 GSE13547_CTRL_VS_ANTI_… 0.00176 0.00855 -0.721 -3.44            0   174 down 
##  2 GSE15750_DAY6_VS_DAY10… 0.00177 0.00855 -0.712 -3.41            0   181 down 
##  3 GSE14415_NATURAL_TREG_… 0.00176 0.00855 -0.689 -3.28            0   173 down 
##  4 GSE15750_DAY6_VS_DAY10… 0.00177 0.00855 -0.676 -3.26            0   188 down 
##  5 GSE39556_CD8A_DC_VS_NK… 0.00176 0.00855 -0.667 -3.22            0   185 down 
##  6 GSE24634_TEFF_VS_TCONV… 0.00176 0.00855 -0.660 -3.19            0   192 down 
##  7 GSE25088_WT_VS_STAT6_K… 0.00178 0.00855 -0.664 -3.17            0   176 down 
##  8 GSE30962_PRIMARY_VS_SE… 0.00177 0.00855 -0.653 -3.14            0   184 down 
##  9 GSE36476_CTRL_VS_TSST_… 0.00178 0.00855 -0.655 -3.14            0   183 down 
## 10 GSE12845_IGD_POS_BLOOD… 0.00176 0.00855 -0.646 -3.12            0   186 down 
## # … with 20 more rows
fgsea_plot(fgsea_res$res_msg, pathways_title='All signatures', condition_name= cond_nm)

## # A tibble: 30 x 8
##    pathway                     pval   padj     ES   NES nMoreExtreme  size state
##    <chr>                      <dbl>  <dbl>  <dbl> <dbl>        <dbl> <int> <chr>
##  1 SOTIRIOU_BREAST_CANCER_… 0.00187 0.0152 -0.818 -3.78            0   138 down 
##  2 ROSTY_CERVICAL_CANCER_P… 0.00185 0.0152 -0.820 -3.75            0   130 down 
##  3 SHEDDEN_LUNG_CANCER_POO… 0.00172 0.0152 -0.679 -3.60            0   426 down 
##  4 WHITEFORD_PEDIATRIC_CAN… 0.00181 0.0152 -0.785 -3.50            0   108 down 
##  5 GNF2_CCNA2               0.00188 0.0152 -0.867 -3.50            0    65 down 
##  6 DUTERTRE_ESTRADIOL_RESP… 0.00182 0.0152 -0.680 -3.49            0   305 down 
##  7 KONG_E2F3_TARGETS        0.00180 0.0152 -0.810 -3.49            0    88 down 
##  8 GNF2_CCNB2               0.00187 0.0152 -0.886 -3.48            0    54 down 
##  9 MODULE_54                0.00186 0.0152 -0.697 -3.48            0   241 down 
## 10 GSE13547_CTRL_VS_ANTI_I… 0.00180 0.0152 -0.721 -3.45            0   174 down 
## # … with 20 more rows